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ORIGINAL ARTICLE 

A Mechanistic, Model-Based Approach to Safety 
Assessment in Clinical Development 

J Lippert 1 , M Brosch 2 , 0 von Kampen 2 , M Meyer 1 , H.-U Siegmund 1 , C Schafmayer 3 , T Becker 3 , B Laffert 4 , L Gorlitz 1 , S Schreiber 2 , 
PJ Neuvonen 56 , M Niemi 56 , J Hampe 2 and L Kuepfer 1 

Assessing the safety of pharmacotherapies is a primary goal of clinical trials in drug development. The low frequency of 
relevant side effects, however, often poses a significant challenge for risk assessment. Methodologies allowing robust 
extrapolation of safety statistics based on preclinical data and information from clinical trials with limited numbers of patients 
are hence needed to further improve safety and efficacy in the drug development process. Here, we present a generic systems 
pharmacology approach integrating prior physiological and pharmacological knowledge, preclinical data, and clinical trial 
results, which allows predicting adverse event rates related to drug exposure. Possible fields of application involve high-risk 
populations, novel drug candidates, and different dosing scenarios. As an example, the approach is applied to simvastatin 
and pravastatin and the prediction of myopathy rates in a population with a genotype leading to a significantly increased 
myopathy risk. 
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Investments in pharmaceutical research and development 
have been on the rise continuously for decades. At the same 
time, the number of annual drug approvals shows a tendency 
to decrease to a level insufficient to ensure a sustainable eco- 
nomic basis for pharmaceutical innovator companies. Recent 
estimates for the total cumulated costs of a newly approved 
drug reach values in the order of US$2 billion. 1 Today, the 
typical duration of the clinical development phase of a drug 
is 8 years and accounts for more than 50% of the total prod- 
uct costs on research and development. 1 The key driver of 
development costs is project attrition due to insufficient clinical 
results. Although the lack of efficacy dominates the statistics of 
attrition causes, a stunning 20% of clinical drug failures result 
from safety issues 2 and numerous approved drugs had to be 
withdrawn from the market due to safety issues arising after 
marketing authorization. 3 The low frequency of relevant side 
effects, however, often poses a significant challenge for risk 
assessment. An improved understanding and predictability of 
the occurrence of adverse events would therefore enable sig- 
nificant improvements in drug development and patient safety. 

In this study, we present a generic approach to exploit drug 
safety-related information routinely generated during preclinical 
and clinical drug development in a quantitative way (Figure 1). 
Information used may also originate from classical in silico toxi- 
cology approaches such as quantitative structure-activity rela- 
tionship, descriptor-, and rule-based models. 4-6 The initial steps 
of our proceeding are state of the art in environmental toxicology 
modeling. 7 In addition, our approach accounts for the existing 
knowledge on genetic risk factors of excess drug exposure or 
drug response as, for example, identified in association studies 8 
and clinical trial information. To identify and quantify potential 



safety issues in high-risk patient populations, we propose to 
integrate available knowledge and prior information together 
with experimental data originating from the drug development 
project in one unified computational model representation, the 
so-called physiology-based pharmacokinetic (PBPK) models 
(Supplementary Figure S1 online). From the establishment 
of the fundamental concepts, 9 PBPK modeling has increas- 
ingly gained acceptance and is nowadays well established in 
environmental toxicology and risk assessment 10-12 as well as 
in drug development. 13-18 PBPK models integrate prior ana- 
tomical and physiological information ranging from the whole 
body level (e.g., organ volumes, blood flow rates, tissue com- 
position) 13 ' 1920 to relative tissue-specific gene expression data 
for relevant metabolic enzymes (e.g., cytochrome P450 3A4) 
and transporters (e.g., solute carrier organic anion transporter 
family member 1B1 (SLC01B1)) at the cellular scale. 21 Like- 
wise, substance-specific absorption, distribution, metabolism, 
and excretion (ADME) properties such as molecular weight, 
lipophilicity, protein binding, or metabolic stability are consid- 
ered for parameterization of the underlying generic distribution 
models. 13 ' 1920 

The experimental data used in our case study range from 
the preclinical characterization of enzymes and transporter 
proteins metabolizing and transporting the drug to safety event 
rates in clinical trials (Figure 1). As PBPK models explicitly 
represent ADME genes, they are also capable of describing 
patient groups and individuals with specific genotypes and 
their corresponding pharmacokinetic phenotype. Another 
unique feature of PBPK models is the explicit representation 
of tissue, thereby, enabling prediction of concentration-time 
profiles based on plasma PK data. The systematic integration 
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Figure 1 Systems pharmacology. A generic systems pharmacology approach exploits drug safety-related information to predict adverse 
event rates. Physiology-based pharmacokinetic models are used to represent large amounts of prior physiological, anatomical, and 
pharmacological information and knowledge. The workflow additionally integrates preclinical data and clinical trial results routinely 
generated drug development. Together with the existing knowledge of genetic risk factors, the approach quantifies potential safety issues 
in high-risk patient populations. Gl, gastro-intestinal; RD, embryonal rhabdomyosarcoma cells; SNP, single-nucleotide polymorphism. 

of prior (pre-)clinical knowledge and information enables as 
such a quantitative prediction of drug exposure in target tis- 
sue for the estimation of toxic effects and the mechanistic 
analysis of adverse drug reactions. 22 
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RESULTS 

The proposed PBPK modeling-based workflow consists of 
six consecutive steps that allow the prediction of adverse 
event rates related to drug exposure (Figure 2, see Methods, 
Supplementary Table S1 online). To demonstrate the feasi- 
bility of our approach, the case of statin-induced myopathy is 
analyzed here. Statins are 3-hydroxy-3-methyl-glutaryl-CoA 
reductase inhibitors that are generally well tolerated and are 
used in the treatment of more than 30 million patients world- 
wide. Although mild cases of myopathy occur in around 1-5% 
of statin-treated patients, only 0.001 % develop rhabdomyolysis, 
i.e., muscle symptoms with a more than 10-fold increase in cre- 
atinine kinase 23 and with the potential for fatal consequences. 
Given the extent of statin use, a relevant group of hundreds of 
patients is therefore at risk for myopathy. In earlier studies, a 
single-nucleotide polymorphism (SNP;c.521T->C, p.Val174Ala) 
in the SLC01B1 encoding for the organic anion transporting 
polypeptide OATP1 B1 has been linked to an increased risk of 
myopathy after simvastatin treatment. 8 ' 2425 For unknown rea- 
sons, the myopathy risk of smaller doses of simvastatin or of cer- 
tain other statins such as pravastatin is significantly lower. 2426 

PBPK model establishment 

In step 1, PBPK models for the clinical pharmacokinetics of 
simvastatin, our reference compound, and pravastatin, our 
assessment candidate, were established in volunteers with 
characterized SLCOW1 genotypes 2728 (see Methods and Sup- 
plementary Materials and Methods online). Simvastatin is a 
lipophilic compound administered as an inactive lactone prod- 
rug. The active acid metabolite is formed by esterase-mediated 



2) Model 
evaluation 




Figure 2 Workflow for model-based safety assessment. (1) Esta- 
blishment of reference PBPK models; (2) model evaluation at 
relevant scales; (3) simulation of virtual populations and model 
evaluation; (4) calculation of TD markers; (5) evaluation of the safety 
risk; (6) prediction of the safety risk in high-risk patient groups or 
the risk for a novel drug candidate as a (a) dose to dose, (b) drug 
to drug, or (c) patient to patient extrapolation. ADME, absorption, 
distribution, metabolism, and excretion; PBPK, physiology-based 
pharmacokinetic; TD, toxicodynamic. 

hydrolysis in the plasma and liver. 29 Both the lactone prodrug 
and the active acid were explicitly represented in the simvas- 
tatin PBPK model. Pravastatin is considerably less lipophilic 
and is partially metabolized by sulfotransferases in the intes- 
tinal epithelium and liver. 30 In addition to OATP1B1, which is 



CPT: Pharmacometrics & Systems Pharmacology 



Model-based safety assessment 

Lippert etal 



represented in the simvastatin model, pravastatin is also trans- 
ported by MRP2 in the intestine, liver and kidney, and by renal 
OAT3. 31 Gene expression data were used as surrogates for 
tissue-specific transporter activity (i.e., relative V max values, see 
Methods). 21 Due to the large amount of prior anatomical and 
physiological information provided in PBPK models, only a lim- 
ited set of parameters had to be adjusted in both models for the 
reference genotype (C.521TT). Here, plasma concentration pro- 
files for oral administration of simvastatin 28 and for intravenous 
and oral administration of pravastatin 273233 were used. In addi- 
tion, urinary excretion data were considered to close the mass 
balance for the compounds. 3435 Both statin PBPK models for 
the reference genotype (C.521TT) are compared with the exper- 
imental data (Figure 3a,b, Supplementary Figures S2-S5 
online, Supplementary Tables S2-S6 online). 

Subsequently, both models were adjusted to the increased 
plasma concentration levels for the homozygous genotype 
harboring the risk allele C.521C by means of a study in 31-32 
healthy Caucasian volunteers with a characterized SLC01B1 
c.521T^C genotype. 2728 We assumed that decreasing 
OATP1B1 transporter activity alone should be sufficient to 
describe the pharmacokinetic difference between the TT and 
the CC group. The simvastatin and pravastatin models for the 
C.521CC genotype were therefore independently fitted based 
on the C.521TT models by adjusting OATP1B1 transporter 
activity only. The excellent representation of the clinical data 
obtained with this approach confirms the assumption of the 
dominating role of the transporter for the difference between the 
clinical phenotypes of the two homozygous volunteer groups 
(Figure 3a,b). Most notably, maximum plasma and muscle tis- 
sue concentrations and integral drug plasma exposure increase 
by a factor of 2-4 in volunteers with the CC genotype as com- 
pared with those with TT genotype for both statins. 36 

Model evaluation 

Next, the predictive power of our models for pharmacokinetic 
phenotypes was evaluated. For this purpose, the PK of the 
heterozygous genotype was predicted by taking the arithmetic 
mean of transporter activities in both homozygous groups. The 
simulated plasma concentration curves of simvastatin and 
pravastatin were compared with the experimentally determined 
PK in the heterozygous genotype C.521TC (Figure 3a,b). 
Visual check indicates the accurate representation of the 
experimental data by the models and the predictive power of 
the computational model for pharmacokinetic phenotypes. 

To also validate our model at the molecular scale, the effect 
of the c.521T->C SNP on kinetic transporter efficiency was 
investigated. To take into account the variability of model 
parameters, an ensemble modeling approach 37 was performed 
based on which the ratio of the genotype-specific transporter 
activities (r = V^JV^J was found t0 be r = °- 265 ( Cl: (0-263; 
0.267)) for simvastatin acid and r = 0.11 (CI: (0.105; 0.118)) 
for pravastatin (Figure 3c). Of note, the transporter protein 
abundance is canceled out when the ratios of the genotype- 
specific transporter activities for pravastatin and simvastatin 
are divided by each other. In this case, the resulting quotient 
of ratios is a mere function of the kinetic transport efficiencies 
(see Methods). The ratios for both statins are statistically 
significantly different (randomization test, P < 1.0 x 10~ 4 ) 
suggesting a substrate-specific effect of the transporter 
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Figure 3 Validation of physiology-based pharmacokinetic (PBPK) 
models for simvastatin acid and pravastatin at the organism level 
and at the molecular scale. (a,b) Model-based prediction of phar- 
macokinetic phenotypes. After adjustment of model parameters 
with respect to the homozygous genotype TT (black; experiment: 
triangle (up)), the minor frequent homozygous genotype (CC, light 
gray; experiment: triangle (down)) was simulated by decreasing the 
transporter activity for (a) simvastatin acid 28 and (b) pravastatin 27 
(all simulations are indicated by a solid line). Taking the average of 
the transporter activity of both homozygous genotypes correctly 
predicts the plasma curves of the heterozygous genotype (TC, dark 
gray; triangle diamonds) for both drugs, (c) Model-based prediction 
of the effect of the solute carrier organic anion transporter family 
member 1B1 single-nucleotide polymorphism on the transporter 
activities at the molecular scale. The ratio of transporter activities 
( r = KSS/K^x) is different for simvastatin acid (black) and pravasta- 
tin (white). This PBPK model-based finding (left) could be verified 
with in vitro assays (right) using HEK 293 cells and [ 3 H]-labeled 
pravastatin and simvastatin acid. *P< 1e-4, randomization test. 



genotype on the kinetic transport efficiencies. 

To validate this model-based finding experimentally, in vitro 
experiments with transfected HEK 293 cells and [ 3 H]-labeled 
pravastatin and simvastatin acid were performed to mea- 
sure the genotype-specific transport rate (see Methods). 
The experimental ratio of transport activities for simvastatin 
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(0.62, CI: (0.1 7; 1 . 1 1 )) is two times higher than for pravastatin 
(0.32, CI: (0.10; 0.51)) (Figure 3c, randomization test, P< 1 
x 10 -4 ; Supplementary Figure S6 online, Supplementary 
Table S8 online). The quotient of transport ratios is close to 
the quotient obtained from the modeling of clinical data and 
shows the excellent capability of the model to quantitatively 
relate in vivo to in vitro observations. This also shows that 
the c.521T^C SNP modifies the kinetic transport efficiency 
which, in turn, contributes to the difference between the clini- 
cal phenotypes. Both PBPK models were hence capable of 
predicting the pharmacokinetic phenotypes and the effect of 
the SLCOW1 SNP at the molecular scale and were therefore 
considered as validated at the relevant biological scales. 

Calculation of a toxicodynamic marker in virtual patient 
populations 

In step 3, Monte-Carlo simulations of genotype-specific virtual 
patient populations were performed of 1,000 individuals each 
by varying the anatomical and physiological parameters in 
the mean PBPK models according to prior statistical informa- 
tion 3839 (see Methods and Supplementary Table S7 online). 
The simulated PK of virtual populations for the three geno- 
types was compared with clinical data 2728 (Figure 4a-f). 

In step 4, systemic effects of the c.521T->C polymorphism 
on exposure in the target tissue were integrated with in vitro 
toxicity to derive a toxicodynamic (TD) marker. Simulations of 
the maximum concentrations for the CC genotype in muscle 
tissue (0.0035 and 0.1 umol/l for simvastatin acid and pravas- 
tatin, respectively) show a 29-fold higher exposure for pravas- 
tatin as compared with simvastatin. By dividing these predicted 
tissue concentrations by half maximal inhibitory concentrations 
values determined with embryonal rhabdomyosarcoma cells 
(3.99 umol/l and 4890 umol/l for simvastatin acid and pravas- 
tatin, respectively), 40 an in vivo marker for statin toxicity was 
calculated. Notably, this TD marker is considerably higher for 
simvastatin (8.8 x 10 -4 ) than for pravastatin (2.1 x 10 -5 ). 

Prediction of myopathy incidence rates in different 
patient populations 



In step 5, cumulated distributions of the TD marker for both 
drugs (Figure 5a-c) were calculated based on the population 
simulations performed earlier (Figure 4a-f). The population 
simulations show a significantly higher risk after simvastatin 
treatment as compared with pravastatin (in agreement with 
clinical myopathy frequencies). 2426 Clinical risk data were 
used to consecutively predict safety event rates both of high- 
risk patient populations and of the assessment candidates 
simultaneously investigated in the workflow. To translate the 
TD marker to the myopathy incidence rates, a range of the TD 
marker was identified such that the fraction of patients having a 
larger value agrees simultaneously with the clinically observed 
first-year incidence rates of TT and TC patients receiving 80 
mg simvastatin in the SEARCH study. 8 By using the binomial 
confidence intervals for the clinical rates, a threshold corridor 
for our TD marker between 0.002 and 0.0024 was thereby 
derived (Figure 6). This corridor enables the drug-drug and 
patient-patient extrapolation, as predicted incidence rates 
correspond to the fraction of simulated patients having a TD 
marker exceeding the established threshold corridor. 

In the final step, the safety risk during various scenarios of 
statin treatment was evaluated. First, the safety risk in a differ- 
ent, high-risk patient subgroup of CC patients was predicted. 
Taking the intersection of the threshold corridor with the simu- 
lated cumulative distribution of the TD marker for the CC geno- 
type, the model predicts an incidence rate between 22.2 and 
29.6% in the first treatment year. Of note, an incidence rate 
of 15.0% (CI: (8.8%; 23.1%)) was reported in the SEARCH 
study. We next simulated the myopathy incidence rate for 40 
mg simvastatin in the Heart Protection Study 41 (see Methods, 
Supplementary Figure S7 online). As the patient population of 
this study was not genotyped for the SLCOW1 SNP, a weighted 
sum of incidence rates was calculated using empirical allele fre- 
quencies of the two SLCOW1 variants. 42 The predicted overall 
myopathy incidence rate between 0.04 and 0.1 1% also agrees 
well with 0.097% found in the study. 41 Finally, the safety risk for 
the putative assessment candidate pravastatin was predicted 
(see Supplementary Figure S8 online). We found no intersec- 
tion of the threshold corridor identified from clinical data of the 
SEARCH study 8 with the simulated cumulative distribution of 
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Figure 5 A toxicodynamic (TD) marker for statin toxicity, (a) The TD marker for pravastatin (dashed line) and simvastatin acid (solid line) is 
plotted as a function of transporter activity relative to the mean transport activity of the TT genotype. Mean activities for all genotypes are 
indicated by circles (TT: black, TC: dark gray, CC: light gray) and the extrapolated intermediate values are indicated by a solid line. (b,c) 
Genotype-specific distributions of the TD marker in population PK studies (TT: black bar, TC: dark gray bar, CC: light gray bar) for 1,000 
virtual patients (cumulated distribution, black line) on treatment with (b) simvastatin acid and (c) pravastatin. The distributions are based 
on the population PK plots in Figure 4. 




Figure 6 Prediction of clinical incidence rates. Simulated cumulated 
distributions of the toxicodynamic (TD) marker are shown as solid 
lines (TT: black, TC: dark gray, CC: light gray). First-year myopathy 
incidence rates of TT and TC patients (vertical ranges, right y axis) 
receiving 80 mg simvastatin in the SEARCH study 8 combined with 
the simulated cumulative distribution of the TD marker are used to 
determine a threshold corridor (vertical light gray bar) for the simulated 
cumulative distribution of the TD marker in different subgroups of 
patients (dashed horizontal lines). The corresponding intersections 
of the threshold corridor and the simulated cumulated distributions 
of the CC genotype are used to forecast the clinical incidence rate 
found in the SEARCH study 8 (horizontal light gray bar). 

the TD marker for pravastatin and predict a vanishing geno- 
type-associated risk for pravastatin. 

DISCUSSION 

Here, we present a generic systems pharmacology frame- 
work approach using PBPK modeling and prior physiologi- 
cal and pharmacological knowledge that may be used to 
extrapolate from clinical incidence rates in dominant patient 
populations to rare frequency events, i.e. patient to patient, 
dose to dose, or drug to drug. In previous publications, it 



was already demonstrated that PBPK models are capable 
of accurately representing the pharmacokinetic phenotypes 
resulting from differences in polymorphic ADME genes linked 
to codeine toxicity. 1443 In these studies, however, the analy- 
sis was restricted to systemic exposure to the drug and its 
metabolites. The development of a generic workflow translat- 
ing genotypes to clinical incident rates in a mechanistic and 
quantitative manner is a fundamental extension. 

The use of our workflow is demonstrated by quantitative 
extrapolation of myopathy rates for the high-risk SLCOW1 
genotype under simvastatin treatment based on incidence 
rates in low-risk reference populations. Likewise, we pre- 
dicted the preferable safety profile for a smaller dose of 
simvastatin. Extrapolations based on both the models are in 
general quantitative agreement (contained within CIs) with 
findings from the SEARCH study 8 and the Heart Protection 
Study. 41 With regard to drug to drug extrapolation, clinical 
simvastatin data were used to predict pravastatin myopathy 
rates. Here, a superior safety profile was estimated. Although 
this is a rather qualitative finding, it is again in line with clinical 
trial data for pravastatin that show no genotype-associated 
risk. 24 Recent pharmacovigilance data of the FDA show 212 
co-occurrences of pravastatin therapy and myopathy, 44 how- 
ever. A baseline myopathy risk resulting from other factors 
such as comedication or renal impairment may explain this 
apparent difference to our prediction, as the model used in 
this study does not account for risk factors other than the 
SLCOW1 SNR This is a conscious restriction of our study 
as the reported trial data do not document additional risk 
factors or the baseline myopathy risk in the general popula- 
tion. In principle, PBPK models are able to represent realis- 
tic combinations of differences in renal function, drug-drug 
interactions, and genotypes as shown earlier. 43 The far lower 
incidence rate of such secondary factors is obviously masked 
by the dominant number of SCLO 1B 7-related events in the 
case of simvastatin. Notably, 2,278 cases are reported for 
simvastatin 44 confirming the predicted superior safety of 
pravastatin in clinical practice, too. We conclude that our 
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approach contributes significantly to a quantitative, mecha- 
nistic understanding of the dramatically different myopathy 
rates between simvastatin and pravastatin. 

A key step of our workflow is the model validation in step 
2, which secures validity of the consecutive toxicity-related 
predictions. Although a formal validation of PKPK models is 
ultimately not possible, the successful prediction of a suf- 
ficiently diverse set of experiments that has not been used 
for model establishment itself generates confidence in the 
quality of the models. In our study, targeted experimental 
data at the molecular and the organism scale were used to 
carefully verify the predictive accuracy of both models. In 
this regard, it is clearly desirable to use as much experimen- 
tal information as possible to support model identification 
and to further strengthen the confidence in the models. It 
should be noted that our model predictions may be assumed 
to be robust as parameter uncertainty was explicitly taken 
into account by ensemble modeling 37 during model estab- 
lishment. In general, the availability of specific experimental 
data and knowledge is a minimal requirement for the appli- 
cability of our approach and also a potential limitation. The 
minimal set consists of ADME properties of a drug, plasma 
PK data, and knowledge about the existence, and the target 
tissue of a toxic effect. 

Multiple in silico methodologies are established tools for 
the prediction of toxicological mechanisms and TD param- 
eters, such as half maximal inhibitory concentrations, in 
drug development. 4 As in the environmental toxicology appli- 
cations of PKPB modeling, 7 ' 1112 our study uses a detailed, 
mechanistic PBPK model to derive toxicity-related predic- 
tions. We use the PBPK model for a unified representation of 
knowledge generated during the drug development process 
that allows, for example, a model-based extrapolation of in 
vitro data to an in vivo context. In our study, this was done, 
for example, for the calculation of the TD marker from in vitro 
half maximal inhibitory concentrations values and the simu- 
lated drug exposure in the target tissue. Going beyond envi- 
ronmental applications, we systematize the steps required to 
integrate clinical information in the process of model estab- 
lishment and validation and the application to predictions for 
specific populations. 

Although investments into the development of novel drug 
candidates have continuously increased for decades, the num- 
ber of annual drug approvals showed a tendency to decline. 
Even worse, adverse events that represent a considerable 
risk for patients at each stage of clinical development still con- 
tribute significantly to the overall clinical attrition rate. Model- 
based approaches that represent all available knowledge and 
information about a drug candidate and foster a mechanistic 
understanding of the pharmacological mechanisms promise to 
improve the success rates in clinical development by enabling 
rational patient stratification and individualized dose selection. 
Using the workflow presented here may help to identify and 
avoid potentially life-threatening side effects. The approach 
presented therefore provides a valuable tool for rational drug 
development. If applied consistently, it may help to significantly 
improve patient safety and support rational decision making at 
crucial stages of clinical development. 



METHODS 

PBPK. PBPK models describe the mechanisms underlying 
the ADME of a substance within the body at in-depth levels 
of detail 1315 16 19 21 384546 (Supplementary Figure S1 online). 
To support model building, PBPK models are based on a 
large degree of prior information regarding an organism's 
anatomy and physiology (see Supplementary Materials 
and Methods online). Such prior physiological information 
may range from the organism level (organ volumes, blood 
flow rate, etc.) to expression profiles of relevant enzymes 
and transporter proteins. Most notably, nearly all model 
parameters are either taken from large-scale collections of 
anatomical and physiological data or are calculated from 
a small set of drug-dependent properties. 13 ' 14 ' 19 " 21 ' 3839 On 
the basis of substance-specific physicochemical surrogate 
parameters such as lipophilicity or plasma protein binding 
(see Supplementary Table S2 online), generic prediction 
models are simultaneously parameterized all over the body 
and can be used to describe drug-concentration profiles 
in various organs and tissues. The basic structure of the 
distribution models describe passive processes based on 
blood flow and diffusion. PBPK models can be extended by 
introducing additional mechanisms such as active trans- 
porters or enzyme-catalyzed clearance processes in vari- 
ous organs. 

The PBPK models for simvastatin and pravastatin were built 
using the software platform consisting of PK-Sim Version 4.2, 
MoBi Version 2.2, and the MoBi Toolbox for Matlab Version 
2.0 (see Supplementary Materials and Methods online). All 
three tools are available for noncommercial academic use free 
of charge (Bayer Technology Services GmbH, Leverkusen, 
Germany). 13 ' 1938 PBPK models with PK-Sim were exported to 
and modified in the MoBi software (Bayer Technology Services 
GmbH, Leverkusen, Germany). All parameter optimizations 
and batch mode simulations were carried out using Matlab 
(version 7; MathWorks, Natick, MA) (see Supplementary 
Materials and Methods online). Both the PBPK models 
can be implemented by using the parameter values given 
in Supplementary Tables S2-S5 online, which provide the 
necessary information. The models are also provided in the 
Supplementary Materials and Methods online (MoBi-files). 

General workflow. The proposed PBPK modeling-based 
workflow consists of six consecutive steps that allow 
predicting adverse event rates related to drug exposure 
(Figure 2, Supplementary Table S1 online). In the first 
step, PBPK models for a frequent reference PK phenotype 
(usually represented by the population mean PK of the most 
frequent patient group) are established for all drugs and 
drug candidates under investigation. In the second step, the 
PBPK models are evaluated to secure validity of the con- 
secutive toxicity-related predictions. In the third step, the PK 
of the reference drug and further assessment candidates is 
simulated for a large virtual reference patient simulation and 
additional virtual patient populations of interest. For each 
simulated virtual individual, plasma PK and target tissue 
concentration-time profiles are calculated. In step 4, the lat- 
ter are integrated with in vitro toxicity information to calculate 
TD markers for all individuals. In step 5, clinical risk data for 
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the dominant reference population and the reference drug 
are introduced as the anchor point for the prediction of safety 
risks of high-risk populations and of the other drugs included 
in the workflow. The final step 6 predicts the safety risk in 
high-risk patient groups or the risk for a novel drug candidate 
as a patient to patient, dose to dose, or drug to drug extrapo- 
lation (Figure 2, Supplementary Table S1 online). 

Ratio of transporter activities. The transport activity is the 
product of the kinetic transport efficiency k cat and the overall 
concentration of OATP1B1 (E 0 ) 

V = k F 

"max A cat *-0 

Both parameters are potentially specific for each allele. Divid- 
ing the ratios 



for both pravastatin and simvastatin reveals that protein abun- 
dance cancels out and that the quotient of ratios is merely 
dependent on the kinetic transport efficiency k \ 
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Transport assay. Pravastatin and simvastatin lactone were 
purchased from Sigma-Aldrich (St. Louis, MO). Simvastatin 
lactone was converted into the acid using a method previously 
described. 47 [ 3 H]-labeled pravastatin and simvastatin acid were 
purchased from Hartmann Analytic (Braunschweig, Germany). 
After transfection for 48 h, the cells were washed with pre- 
warmed phosphate-buffered saline and incubated in fetal 
calf serum (FCS)-free Dulbecco's modified Eagle's medium 
(DCEM). After 2 h, the culture medium was replaced with FCS- 
free DMEM supplemented each with pravastatin and activated 
simvastatin in a concentration of 20 and 10 umol/l, respec- 
tively, each spiked with 0.5 uCi of the respective [ 3 H]-labeled 
substrate. These concentrations are lower than those at which 
the half maximal uptakes occur. 48 After 30-min incubation in 
5% C0 2 at 37 °C, the cells were washed with ice cold DMEM 
supplemented with 5% FCS twice and additionally three times 
with ice cold FCS-free DMEM. Subsequently, cells were lysed 
with 800 ul 0.2% sodium dodecyl sulfate in phosphate-buffered 
saline. Uptake of labeled substrate was measured by liquid 
scintillation counting. Each sample was measured at least in 
triplicates and normalized to total protein concentration (see 
Supplementary Materials and Methods online). 

Calculation of confidence intervals. We have reported 95% 
confidence intervals that have been calculated by ran- 
domization. Calculation was performed by estimating the 
sampling distribution using permutation of inputs (cf. per- 
mutation tests). This gives the confidence intervals under 
the null hypothesis and that the parameters are indepen- 
dent of the outcome. The experimentally measured data are 
assumed to be normally distributed around the estimated 
mean with estimated SD. The confidence interval for the 



mean is derived using standard formulas. All statistics for 
derived quantities are calculated by generating uniformly 
distributed values from the respective confidence intervals 
and by evaluating the quantity. The lower and upper limits of 
the reported confidence intervals are the 2.5%-quantile and 
the 97.5%-quantile of this computationally generated prob- 
ability distribution. 

Prediction of clinical incidence rates (Heart Protection 
Study), 40 mg dose of simvastatin. The threshold on the TD 
marker identified for patients with TT and TC genotype in the 
SEARCH study was then applied to risk score data calcu- 
lated in simvastatin PBPK models following 40 mg dosings of 
simvastatin. Here, an incidence rate is calculated per geno- 
type, which was then weighted by the measured genotype 
frequency in the Heart Protection Study 41 (Supplementary 
Figure S7 online). 
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Study Highlights 

WHAT ISTHE CURRENT KNOWLEDGE ON 
THETOPIC? 

Asingle-nucleotide polymorphism in SLCOW1 has been 
linked to an increased risk of myopathy after simvastatin 
treatment. 

WHAT QUESTION DIDTHIS STUDY ADDRESS? 

It is unknown, however, why the myopathy risk of certain 
other statins such as pravastatin is significantly lower. 

WHATTHIS STUDY ADDSTO OUR KNOWLEDGE 

Using an integrative systems pharmacology approach 
of modeling and quantitative experiments, we first dem- 
onstrate that the specific SLCOW1 SNP directly affects 
kinetic transport efficiency at the molecular scale. The 
model was then used to calculate genotype-specific 
TD markers based on in vitro measurements for statin 
cytotoxicity. The model could thereby predict clinical inci- 
dence rates of myopathy in high-risk subpopulation for 
novel drug candidates and different dosing scenarios. 

HOWTHIS MIGHT CHANGE CLINICAL 
PHARMACOLOGY ANDTHERAPEUTICS 

In this study, we show how unified computational models 
integrating prior knowledge together with specific experi- 
mental data can be used to identify and quantify safety 
risks prior to large clinical trials and may, hence, signifi- 
cantly support rational drug design. 
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